EEA 2020 - TP2 - Darío Hruszecki

# Carga de librerías
library(tidyverse)

## ── Attaching packages ─────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────── tidyverse 1.3.0 ──

## ✓ ggplot2 3.3.2     ✓ purrr   0.3.4
## ✓ tibble  3.0.3     ✓ dplyr   1.0.0
## ✓ tidyr   1.1.0     ✓ stringr 1.4.0
## ✓ readr   1.3.1     ✓ forcats 0.5.0

## ── Conflicts ────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────── tidyverse_conflicts() ──
## x dplyr::filter() masks stats::filter()
## x dplyr::lag()    masks stats::lag()

library(tidymodels)

## ── Attaching packages ────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────── tidymodels 0.1.1 ──

## ✓ broom     0.7.0      ✓ recipes   0.1.13
## ✓ dials     0.0.9      ✓ rsample   0.0.8 
## ✓ infer     0.5.3      ✓ tune      0.1.1 
## ✓ modeldata 0.1.0      ✓ workflows 0.2.1 
## ✓ parsnip   0.1.4      ✓ yardstick 0.0.7

## ── Conflicts ───────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────── tidymodels_conflicts() ──
## x scales::discard() masks purrr::discard()
## x dplyr::filter()   masks stats::filter()
## x recipes::fixed()  masks stringr::fixed()
## x dplyr::lag()      masks stats::lag()
## x yardstick::spec() masks readr::spec()
## x recipes::step()   masks stats::step()

library(epiDisplay)

## Loading required package: foreign

## Loading required package: survival

## Loading required package: MASS

## 
## Attaching package: 'MASS'

## The following object is masked from 'package:dplyr':
## 
##     select

## Loading required package: nnet

## 
## Attaching package: 'epiDisplay'

## The following object is masked from 'package:yardstick':
## 
##     kap

## The following object is masked from 'package:scales':
## 
##     alpha

## The following object is masked from 'package:ggplot2':
## 
##     alpha

library(GGally)

## Registered S3 method overwritten by 'GGally':
##   method from   
##   +.gg   ggplot2

library(corrr)
library(knitr)
library(ggrepel)
library(kableExtra)

## 
## Attaching package: 'kableExtra'

## The following object is masked from 'package:dplyr':
## 
##     group_rows

library(here)

## here() starts at /home/daro/cdatos/TPs/EEA/TP2

EEA 2020 - TP2 - Darío Hruszecki

DATOS

El dataset de precios de inmuebles proviene de Properati. El mismo ya fue filtrado por los docentes y a su vez se encuentra particionado en subconjuntos de training y test.

Carga partición de entrenamiento

train_ds <- read_csv(here("/ds/ar_properties_train.csv"))

## Parsed with column specification:
## cols(
##   id = col_character(),
##   l3 = col_character(),
##   rooms = col_double(),
##   bathrooms = col_double(),
##   surface_total = col_double(),
##   surface_covered = col_double(),
##   price = col_double(),
##   property_type = col_character()
## )

Describimos el dataset

Previsualizamos los datos

head(train_ds)

Vemos las dimensiones del dataset

dim_desc(train_ds)

## [1] "[32,132 x 8]"

Vemos como se distribuyen los valores de las variables

summary(train_ds)

##       id                 l3                rooms         bathrooms    
##  Length:32132       Length:32132       Min.   :1.000   Min.   :1.000  
##  Class :character   Class :character   1st Qu.:2.000   1st Qu.:1.000  
##  Mode  :character   Mode  :character   Median :3.000   Median :1.000  
##                                        Mean   :2.722   Mean   :1.427  
##                                        3rd Qu.:3.000   3rd Qu.:2.000  
##                                        Max.   :8.000   Max.   :5.000  
##  surface_total surface_covered      price        property_type     
##  Min.   : 28   Min.   : 26.00   Min.   : 69500   Length:32132      
##  1st Qu.: 46   1st Qu.: 42.00   1st Qu.:120000   Class :character  
##  Median : 65   Median : 58.00   Median :169000   Mode  :character  
##  Mean   : 80   Mean   : 70.61   Mean   :214991                     
##  3rd Qu.: 99   3rd Qu.: 85.00   3rd Qu.:260000                     
##  Max.   :320   Max.   :265.00   Max.   :950000

Valores únicos para la variable alfanumérica l3

unique(train_ds$l3)

##  [1] "Almagro"              "Palermo"              "Caballito"           
##  [4] "Villa Crespo"         "Villa Urquiza"        "Barracas"            
##  [7] "Nuñez"                "Belgrano"             "Floresta"            
## [10] "Recoleta"             "Colegiales"           "Parque Chas"         
## [13] "Barrio Norte"         "Villa Devoto"         "Villa Ortuzar"       
## [16] "Velez Sarsfield"      "Parque Chacabuco"     "Villa Pueyrredón"    
## [19] "Villa Real"           "Once"                 "Flores"              
## [22] "San Telmo"            "Las Cañitas"          "Villa Santa Rita"    
## [25] "Villa del Parque"     "Retiro"               "Parque Centenario"   
## [28] "Chacarita"            "Abasto"               "San Cristobal"       
## [31] "Liniers"              "Balvanera"            "Boedo"               
## [34] "Villa General Mitre"  "Saavedra"             "Parque Patricios"    
## [37] "Boca"                 "Paternal"             "Monserrat"           
## [40] "San Nicolás"          "Parque Avellaneda"    "Centro / Microcentro"
## [43] "Puerto Madero"        "Congreso"             "Monte Castro"        
## [46] "Mataderos"            "Coghlan"              "Versalles"           
## [49] "Villa Lugano"         "Catalinas"            "Villa Luro"          
## [52] "Constitución"         "Pompeya"              "Tribunales"          
## [55] "Agronomía"            "Villa Soldati"        "Villa Riachuelo"

Valores únicos para la variable alfanumérica property_type

unique(train_ds$property_type)

## [1] "Departamento" "PH"           "Casa"

Verificamos si hay valores NA en las variables

apply(train_ds, 2, function(x) any(is.na(x)))

##              id              l3           rooms       bathrooms   surface_total 
##           FALSE           FALSE           FALSE           FALSE           FALSE 
## surface_covered           price   property_type 
##           FALSE           FALSE           FALSE

Calculamos la distribución de frecuencias para los valores de property_type

tab1(train_ds$property_type, sort.group = "decreasing", cum.percent = TRUE)

## train_ds$property_type : 
##              Frequency Percent Cum. percent
## Departamento     28203    87.8         87.8
## PH                3120     9.7         97.5
## Casa               809     2.5        100.0
##   Total          32132   100.0        100.0

Calculamos la distribución de frecuencias para los valores de l3

tab1(train_ds$l3, sort.group = "decreasing", cum.percent = TRUE)

## train_ds$l3 : 
##                      Frequency Percent Cum. percent
## Palermo                   4890    15.2         15.2
## Belgrano                  2736     8.5         23.7
## Almagro                   2486     7.7         31.5
## Caballito                 2393     7.4         38.9
## Villa Crespo              2095     6.5         45.4
## Recoleta                  1824     5.7         51.1
## Barrio Norte              1359     4.2         55.3
## Villa Urquiza             1319     4.1         59.4
## Flores                     916     2.9         62.3
## Nuñez                      887     2.8         65.1
## Balvanera                  883     2.7         67.8
## Villa Devoto               527     1.6         69.4
## Colegiales                 523     1.6         71.1
## Villa del Parque           500     1.6         72.6
## San Telmo                  488     1.5         74.2
## Saavedra                   433     1.3         75.5
## San Cristobal              417     1.3         76.8
## Parque Centenario          404     1.3         78.1
## Floresta                   366     1.1         79.2
## Puerto Madero              357     1.1         80.3
## Boedo                      336     1.0         81.3
## Paternal                   327     1.0         82.4
## San Nicolás                325     1.0         83.4
## Monserrat                  310     1.0         84.3
## Chacarita                  292     0.9         85.3
## Retiro                     291     0.9         86.2
## Villa Pueyrredón           285     0.9         87.0
## Parque Chacabuco           282     0.9         87.9
## Barracas                   280     0.9         88.8
## Villa Luro                 264     0.8         89.6
## Liniers                    260     0.8         90.4
## Mataderos                  259     0.8         91.2
## Once                       255     0.8         92.0
## Congreso                   246     0.8         92.8
## Coghlan                    240     0.7         93.5
## Las Cañitas                187     0.6         94.1
## Abasto                     179     0.6         94.7
## Parque Patricios           171     0.5         95.2
## Monte Castro               168     0.5         95.7
## Constitución               139     0.4         96.2
## Villa Santa Rita           128     0.4         96.6
## Villa Lugano               113     0.4         96.9
## Versalles                  113     0.4         97.3
## Villa Ortuzar              109     0.3         97.6
## Centro / Microcentro       105     0.3         97.9
## Boca                       104     0.3         98.3
## Villa General Mitre         99     0.3         98.6
## Parque Chas                 78     0.2         98.8
## Parque Avellaneda           78     0.2         99.0
## Pompeya                     65     0.2         99.2
## Velez Sarsfield             56     0.2         99.4
## Villa Real                  55     0.2         99.6
## Tribunales                  54     0.2         99.8
## Agronomía                   51     0.2         99.9
## Villa Riachuelo             13     0.0        100.0
## Villa Soldati                9     0.0        100.0
## Catalinas                    3     0.0        100.0
##   Total                  32132   100.0        100.0

Observaciones sobre los datos

El dataset contiene 8 variables y 32132 registros
De las 8 variables, 5 son núnericas y 3 son alfanuméricas
De las variables alfanuméricas tenemos el barrio (l3) y el tipo de propiedad (property_type). La tercera es un UUID para la propiedad (id)
En base a los valores únicos encontrados para l3 podemos inferir que las propiedades pertenecen a CABA en su totalidad.
De las númericas, además del precio (price) tenemos variables que describen la propiedad: cantidad de habitaciones (rooms), cantidad de baños (bathrooms), supervicie total (surface_total) y superficie cubierta (surface_covered)
No se observan valores NA en ninguna de las columnas
Aproximadamente un 88% de las propiedades son de tipo Departamento
Aproximadamente un 60% de las propiedades se encuentran en los neighborhood de Palermo, Belgrano, Almagro, Caballito, Villa Crespo, Recoleta, Barrio Norte y Villa Urquiza

Modelo de Regresion Lineal Multiple

Generación del Modelo con todas las covariables (lm_all)

lm_all <- lm(price ~ rooms + bathrooms + surface_total + surface_covered + property_type + l3, data = train_ds)

tidy_lm_all <- tidy(lm_all, conf.int = TRUE)
tidy_lm_all

Análisis del Modelo con todas las covariables

Significado de los coeficientes estimados

(intercept) (-109143.7839): es una inferencia teórica del modelo, no representa la realidad, sería una propiedad sin superficie, habitaciones, barrio, etc.
rooms (-4359.7289): Al ser negativo implica una quita promedio de ese valor en el precio de la propiedad al aumentar en 1 la cantidad de ambientes.
bathrooms (34367.2123): Indica en cuanto aumenta en promedio el valor de la propiedad por agregar un baño.
surface_total (890.5585): Indica en cuanto aumenta en promedio el valor de la propiedad por aumentas en 1 m2 la superficie total.
surface_covered (1497.9310): Indica en cuanto aumenta en promedio el valor de la propiedad por aumentas en 1 m2 la superficie cubierta.
property_typeDepartamento (91485.9410): El valor de β (91485.9410) indica cuanto aumenta la función de respuesta para un departamento en comparación de una casa.
property_typePH (46220.0449): El valor de β (46220.0449) indica cuanto aumenta la función de respuesta para un PH en comparación de una casa.
l3{Barrio} : Son 56 coeficientes que indican en cuanto aumenta o se disminuye en promedio el precio de una propiedad dependiendo del barrio en donde esta ubicada.

Observación: para cada coeficiente se asume que se mantienen constantes el restos de las covariables.

¿Qué observamos al respecto de la significatividad de las variables dummy?

Analizando los p-valores de las variables dummies asociadas a property_type estas resultan estadísticamente significativas. Por otro lado, las variables dummies asociados a la variable l3 se reparten entre algunas significativas y otras no.

Aplicamos Test F

tidy(anova(lm_all))

Analizando la significatividad global de las variables l3 y property_type y observando que tienen un p-valor menor a 0.5 se puede concluir que ambas son signigicativas.

Evaluación del Modelo

Resúmen global (glance) del modelo lm_all

glance(lm_all)

Gráfico de coeficientes del modelo lm_all

ggplot(tidy_lm_all, 
  aes(estimate, term, xmin = conf.low, xmax = conf.high, height = 1)) +
  geom_point(color = "blue", size=3) +
  geom_vline(xintercept = 0, lty = 2, color = "black") + 
  geom_errorbarh(color = "red", size=1) + 
  theme_bw() +
  labs(y = "Coeficientes β", x = "Estimación")

El r.squared en el resúmen global representa el R cuadrado o coeficiente de determinación que refleja la bondad del ajuste del modelo a la variable precio. Al ser un valor cercano a 1 nos está diciendo que el modelo lm_all ajusta muy bien dicha variable.

Por otro lado, el gráfico nos permite visualizar con mayor claridad los intervalos de confiaza de de los coeficientes de las variables dummies de l3. Muchas de ellas incluyen el cero, ergo estas variables podrían considerarse no significativas.

Generación de un modelo lineal sin la covariable l3

lm_wo_l3 <- lm(price ~ rooms + bathrooms + surface_total + surface_covered + property_type, data = train_ds)

tidy_lm_wo_l3 <- tidy(lm_wo_l3, conf.int = TRUE)
tidy_lm_wo_l3

Aplicamos Test F

tidy(
  anova(lm_wo_l3)
)

Al quitar la variable l3 la cantidad de coeficientes se reduce considerablemente, no necesitamos graficar para indentificar que ningúno de los intervalos de confianza incluye el 0. Por otro lado todos los p-values son menores a 0.5 ergo las variables se pueden considerar significativas, esto tambien se ve confirmado por el resultado del Test F.

Determinamos cuál modelo explica mejor la variabilidad del precio (lm_all o lm_wo_l3)

models = list(lm_all = lm_all, lm_wo_l3 = lm_wo_l3)

purrr::map_df(models, broom::glance, .id = "model")

El R cuadrado de lm_all es superior al del modelo lm_wo_l3 por lo tanto el modelo que contiene todas las covariables es el que mejor explica la variable la variabilidad del precio.

Creación de Variables

Creación variable neighborhood

En el análisis anterior se observó que l3 no es un buen agrupador, algunas de sus dummies eran significativas pero muchas otras no. Vamos a proceder a crear una nueva variable neighborhood que agrupe a las propiedad teniendo en cuenta el precio del metro cuadrado.

Es necesario tambien agregar otra variable auxiliar que contenta el precio por metro para luego poder hacer el agrupamiento necesario:

price_m_sq= price_m_sq / surface_total

La valores posibles de neighborhood van a ser:

low_price: price_m_sq <= Q1 (price_m_sq)
medium_price: price_m_sq > Q1 (price_m_sq) & price_m_sq <= Q2 (price_m_sq)
high_price: price_m_sq > Q3 (price_m_sq)

train_ds = train_ds %>% mutate(
    price_m_sq = price / surface_total,
    neighborhood = factor(
      case_when(
        price_m_sq <= quantile(price_m_sq)[2] ~ "low_price",
        price_m_sq > quantile(price_m_sq)[2] & price_m_sq <= quantile(price_m_sq)[4] ~ "medium_price",
        price_m_sq > quantile(price_m_sq)[4] ~ "high_price"
      )
    )
  )

train_ds

Generamos modelo reemplazando l3 or neighborhood

lm_neighborhoods <- lm(price ~ rooms + bathrooms + surface_total + surface_covered + property_type + neighborhood, data = train_ds)

tify_lm_neighborhoods <- tidy(lm_neighborhoods, conf.int = TRUE)
tify_lm_neighborhoods

Determinamos cuál modelo explica mejor la variabilidad del precio (lm_all, lm_wo_l3 o lm_neighborhood)

models = list( lm_all = lm_all, lm_wo_l3 = lm_wo_l3, lm_neighborhoods = lm_neighborhoods )

purrr::map_df(models, broom::glance, .id = "model") %>% arrange(adj.r.squared)

El modelo lm_neighborhoods supera al lm_all en el R cuadrado y por lo tanto tiene mayor explicabilidad sobre la variabilidad del precio de las propiedades.

Generación Modelo Supercicie Descubierta

Dada la fuerte correlación entre las variables surface_total y surface_covered, vamos a generar una nueva variable llamada sup_descubierta que representa la diferencia entre la superficie total y la cubierta.

train_ds = train_ds %>%
  mutate(sup_descubierta = surface_total - surface_covered)

train_ds

Ahora con nuestra nueva variable, vamos a generar un nuevo modelo que incluya la variable sup_descubierta en reemplazo de la variable surface_total.

modelo_sup_descubierta <- lm(
  price ~ rooms + bathrooms + sup_descubierta + surface_covered + property_type + price_m_sq + neighborhood , 
  data = train_ds
) 
tidy_msd <- tidy(modelo_sup_descubierta, conf.int = TRUE)
tidy_msd

Podemos observar como los p-valores de todos nuestros coeficientes resultan significativos en nuestro nuevo modelo.

tidy(
  anova(modelo_sup_descubierta)
)

Y al ejecutar el Test F, validamos el mismo.

Ahora comparamos todos los modelos que generamos.

modelos = list(
  lm_all = lm_all,
  lm_wo_l3 = lm_wo_l3,
  lm_neighborhoods = lm_neighborhoods,
  modelo_sup_descubierta = modelo_sup_descubierta
)

purrr::map_df(modelos, broom::glance, .id = "model") %>%
  arrange(adj.r.squared)

Podemos observar como nuestro nuevo modelo modelo_sup_descubierta es el que mejor explica la variabilidad.

Diagnóstico del Modelo

Vamos a evaluar nuestro modelo lineal con la variable sup_descubierta

plot(modelo_sup_descubierta)

- Residuos vs valores predichos: parece existir una cierta estructura en los datos. La linealidad parece que no se respeta. Vemos la existencía de heterocedasticidad a medida que aumentan las predicciones. Se resaltan 2 observaciones. - Normal QQ plot: los extremos no se ajustan a la distribución teórica, tanto inferior izquierdo como superior derecho. - Residual vs leverage:Existen cuatro puntos con un leverage bastante alto. Y se identifican 2 puntos como posibles outliers. .

Diagnóstico del modelo: El modelo creado no cumple con los supuestos del modelo lineal. No parece ser concluyente respecto de la homocedasticidad, definitivamente existe una falta de normalidad y hay presencia de observaciones de alto leverage

Modelo Log(price)

Vamos a generar un modelo para predecir el logaritmo natural de la variable price \[ log(price) = \beta_0 + \beta_1log(rooms) + \beta_2log(bathrooms) + \beta_3log(surface\_covered) + \beta_4property\_type + \beta_5barrio + \beta_6surface\_patio \]

Para esto vamos a generar nuestros nuevas variables

modelo_log = lm(
  log(price) ~ log(rooms) + log(bathrooms) + log(surface_covered) + property_type + neighborhood + sup_descubierta,
  data = train_ds
)
tidy_ml <- tidy(modelo_log, conf.int = TRUE)
tidy_ml

tidy(
  anova(modelo_log)
)

De nuestro modelo confirmamos que todos los coeficientes resultan significativos.

Ahora comparamos nuestro nuevo modelo_log con nuestro modelo_sup_descubierta.

modelos = list(
  modelo_log = modelo_log,
  modelo_sup_descubierta = modelo_sup_descubierta
)

purrr::map_df(modelos, broom::glance, .id = "model")

Observamos que el mismo explica mejor la variabilidad \(Rˆ2\)

Ahora vamos a analizar el cumplimiento de los supuestos del modelo lineal

plot(modelo_log)

- Residuos vs valores predichos: parece existir una cierta estructura en los datos. La linealidad se respeta. Podemos observar un comportamiento heterocedástico aunque parece ser más uniforme que el del modelo_sup_descubierta. Se resaltan 3 observaciones (observaciones 3636,3709,7346). - Normal QQ plot: los extremos no se ajustan a la distribución teórica, tanto inferior izquierdo como superior derecho. Sin embargo parece que se comportan ligeramente mejor que el modelo_sup_descubierta (parecen estár mas suavizados los extremos). - Residual vs leverage:. se pueden apreciar 2 puntos con un leverage bastante alto y el grafico resalta 3 observaciones.

Diagnóstico del modelo: El modelo creado no cumple con los supuestos del modelo lineal.Sin embargo y de observar los gráficos, parece que se comporta mejor que el modelo_sup_descubierta.

Selección del Modelo

Vamos a generar un nuevo modelo utilizando como variables \(log(rooms)\), \(log(bathrooms)\), \(log(surface_total)\), property_type y neighborhood y otro similar a nuestro modelo_log reemplazando la variable neighborhood por la variable l3.

Generación Nuevos Modelos

Generamos nuestro modelo modelo_log_surface_total que buscará inferir el resultado del \(log(price)\)

modelo_log_surface_total = lm(
  log(price) ~ log(rooms) + log(bathrooms) + log(surface_total) + property_type + neighborhood,
  data = train_ds
)
tidy_mlst <- tidy(modelo_log_surface_total, conf.int = TRUE)
tidy_mlst

En nuestro modelo podemos observar como el coeficiente de \(log(rooms)\) parece no ser significativo (su p-valor = 0.5) y su intervalo de confianza incluye el 0. El resto de nuestros coeficientes parecen ser significativos.

Debido a que neustros coeficientes para la variable \(log(rooms)\) parece que no son significativos, vamos a analizar si la interacción de la variable si lo es para nuestro modelo.

tidy(
  anova(modelo_log_surface_total)
)

Y resulta que si lo es.

Analizamos el cumplimiento de los supuestos de nuestro modelo.

plot(modelo_log_surface_total)

Podemos ver como nuestro modelo se comporta de forma muy similar al modelo_log. Sin embargo hay una diferencia importante al analizar Residuos vs Leverage. Aca podemos observar como hay una especie de corte en el leverage entre los valores de 0.0007 y 0.0013 (aproximadamente).

Generamos nuestro modelo modelo_log_surface_total que buscará inferir el resultado del \(log(price)\)

modelo_log_l3 = lm(
  log(price) ~ log(rooms) + log(bathrooms) + log(surface_covered) + property_type + sup_descubierta + l3 ,
  data = train_ds
)
tidy_mll3 <- tidy(modelo_log_l3, conf.int = TRUE)
tidy_mll3

En este modelo, podemos observar el mismo problema que surgía al incorporar al analisis la variable l3. Sin embargo, el resto de nuestros coeficientes parecen ser estadisticamente significativos.

Debido a que neustros coeficientes para la variable l3 parece que no son significativos, vamos a analizar si la interacción de la variable si lo es para nuestro modelo.

tidy(
  anova(modelo_log_l3)
)

Y resulta que si lo es.

Analizamos el cumplimiento de los supuestos de nuestro modelo.

plot(modelo_log_l3)

- Residuos vs valores predichos: parece existir una cierta estructura en los datos. La linealidad se respeta. Podemos observar un comportamiento heterocedástico (muy similar al modelo_log). Se resaltan 3 observaciones (observaciones 7346,7942, 28022). - Normal QQ plot: los extremos no se ajustan a la distribución teórica, tanto inferior izquierdo como superior derecho. Sin embargo parece haber cierto suavizado en el extremo superior derecho en comparación a otros modelos. - Residual vs leverage:. se pueden apreciar 3 puntos con un leverage bastante alto y el grafico resalta 3 observaciones (3636, 22246,26579).

Diagnóstico del modelo: El modelo creado no cumple con los supuestos del modelo lineal. Parece que se comporta muy similar al modelo_log con una excepción al analizar el grafico de Residual vs leverage.

Evaluación de Modelos Seleccionados

Ahora vamos a elegir 2 de los modelos previamente desarrollados para evaluar y seleccionar, junto a nuestros 2 nuevos modelos, cual es el que mejor permite inferir el precio de una propiedad. Vamos a elegir los modelos que mejor expresan la variabilidad \(Rˆ2\). Estos son el modelo_log y modelo_sup_descubierta.

modelos = list(
  modelo_log = modelo_log,
  modelo_sup_descubierta = modelo_sup_descubierta,
  modelo_log_surface_total = modelo_log_surface_total,
  modelo_log_l3 = modelo_log_l3
)

purrr::map_df(modelos, broom::glance, .id = "model") %>%
  arrange(adj.r.squared)

Podemos observar que el modelo que mejor explica la variabilidad es el modelo modelo_log_surface_total. Observamos el valor del \(Rˆ2\) ajustado.

Predicción

Cargamos el dataset de testing y le aplicamos nuestras transformaciones.

propiedades_test <- read_csv(here("/ds/ar_properties_test.csv"))

## Parsed with column specification:
## cols(
##   id = col_character(),
##   l3 = col_character(),
##   rooms = col_double(),
##   bathrooms = col_double(),
##   surface_total = col_double(),
##   surface_covered = col_double(),
##   price = col_double(),
##   property_type = col_character()
## )

propiedades_test = propiedades_test %>%
  mutate(
    price_m_sq = price / surface_total,
    neighborhood = factor(
      case_when(
        price_m_sq <= quantile(price_m_sq)[2] ~ "low_price",
        price_m_sq > quantile(price_m_sq)[2] & price_m_sq <= quantile(price_m_sq)[4] ~ "medium_price",
        price_m_sq > quantile(price_m_sq)[4] ~ "high_price"
      )
    ),
    sup_descubierta = surface_total - surface_covered
  )

propiedades_test

Evaluamos en Training

Ahora vamos a evaluar nuestros modelos en training y evaluar sus predicciones en testing. Para esto es necesario hacer un análisis distinto para los modelos que utilizan \(log\) que para el modelo_sup_descubierta.

modelos_log = list(
  modelo_log = modelo_log,
  modelo_log_l3 = modelo_log_l3,
  modelo_log_surface_total = modelo_log_surface_total
)

predicciones_training_log = map(.x = modelos_log, .f = augment)

map_dfr(
  .x = predicciones_training_log, 
  .f = rmse, 
  truth = exp(`log(price)`), 
  estimate = exp(.fitted), 
  .id="modelo"
) %>% arrange(.estimate)

eval_train_sd = augment(modelo_sup_descubierta)

rmse(
  data = eval_train_sd,
  truth = price,
  estimate = .fitted
)

Analizando los RMSE de nuestros modelos en training, podemos observar que el modelo con el menor valor es nuestro modelo_sup_descubierta seguido (en orden) por los modelos modelo_log_surface_total, modelo_log, modelo_log_l3.

Evaluamos en Testing

Ahora vamos a evaluar nuestros modelos en testing.

# Aplicamos la función augment a los 4 modelos con el set de testing
predicciones_testing = map(
  .x = modelos_log, 
  .f = augment,
  newdata = propiedades_test
) 

# Obtenemos el RMSE para los 4 modelos
map_dfr(
  .x = predicciones_testing, 
  .f = rmse, 
  truth = price, 
  estimate = exp(.fitted), 
  .id="modelo"
) %>% 
  arrange(.estimate)

pred_log = augment(
  modelo_sup_descubierta, 
  newdata=propiedades_test
) 
pred_log

rmse(
  data = pred_log, 
  truth = price, 
  estimate = .fitted
)

Volvemos a observar como el modelo con el menor valor de RMSE es nuestro modelo_sup_descubierta.

Conclusion

En lineas generales no podemos afirmar que nuestros modelos cumplen los supuestos del modelo lineal.

Sin embargo, podemos utilizar nuestros modelos para predecir los valores del precio de las propiedades. Ahora bien, para elegir el mejor de los modelos propuestos, vamos a fundamentar nuestra decisión en base a la evaluación de las 2 métricas que estamos observando: - \(Rˆ2\) - RMSE

Los modelos que mejor explican la variabilidad medida por \(Rˆ2\) ajustado son el modelo_log, modelo_sup_descubierta y modelo_log_surface_total. Todos estos modelos se encuentran con valores superiores a 0.91.

Sin embargo, nosotros queremos predecir nuevos datos y para ello es importante medir el RMSE para evaluar el error en la predicción. Al analizar esta metrica, observamos un salto importante entre nuestro modelo_sup_descubierta (38974.87) y sus seguidores modelo_log_surface_total (43832.07) y modelo_log (46152.97).

Por ende seleccionaría como mi mejor modelo el modelo modelo_sup_descubierta.

eea2020_tp2_dario_hruszecki

Darío Hruszecki

EEA 2020 - TP2 - Darío Hruszecki

DATOS

Carga partición de entrenamiento

Describimos el dataset

Previsualizamos los datos

Vemos las dimensiones del dataset

Vemos como se distribuyen los valores de las variables

Valores únicos para la variable alfanumérica l3

Valores únicos para la variable alfanumérica property_type

Verificamos si hay valores NA en las variables

Calculamos la distribución de frecuencias para los valores de property_type

Calculamos la distribución de frecuencias para los valores de l3

Observaciones sobre los datos

Modelo de Regresion Lineal Multiple

Generación del Modelo con todas las covariables (lm_all)

Análisis del Modelo con todas las covariables

Significado de los coeficientes estimados

¿Qué observamos al respecto de la significatividad de las variables dummy?

Aplicamos Test F

Evaluación del Modelo

Resúmen global (glance) del modelo lm_all

Gráfico de coeficientes del modelo lm_all

Generación de un modelo lineal sin la covariable l3

Aplicamos Test F

Determinamos cuál modelo explica mejor la variabilidad del precio (lm_all o lm_wo_l3)

Creación de Variables

Creación variable neighborhood

Generamos modelo reemplazando l3 or neighborhood

Determinamos cuál modelo explica mejor la variabilidad del precio (lm_all, lm_wo_l3 o lm_neighborhood)

Generación Modelo Supercicie Descubierta

Diagnóstico del Modelo

Modelo Log(price)

Selección del Modelo

Generación Nuevos Modelos

Evaluación de Modelos Seleccionados

Predicción

Evaluamos en Training

Evaluamos en Testing

Conclusion